// 版权所有2009 Go作者。保留所有权利。
// 此源代码的使用受BSD样式
// 许可证的约束，该许可证可以在许可证文件中找到。

package math

/*
	Floating-point logarithm.
*/

// 下面的原始C代码、长注释和常量
// 来自FreeBSD的/usr/src/lib/msun/src/e_日志。c 
// 并随附此通知。go代码是C.
// 
// ============================================================================================================================================================
// Sun Microsystems，Inc.版权所有。
// 
// 由Sun Microsystems，Inc.的企业SunPro开发。
// 使用、复制、修改和分发本
// 软件的权限是免费授予的，前提是保留本通知
// 的权利。
// /========================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================================。变元约简：求k和f，使
// x=2**k*（1+f），
// ，其中sqrt（2）/2<1+f<sqrt（2）。
// 
// 2。对数（1+f）的近似值。
// 设s=f/（2+f）；基于log（1+f）=log（1+s）-log（1-s）
// /=2s+2/3s**3+2/5s**5+。。。。。，
// /=2s+s*R 
// 我们在[0,0.1716]上使用一个特殊的Reme算法来生成一个14次的多项式来逼近R。这个多项式逼近的最大误差
// 的范围是2**-58.45。也就是说，
// 2 4 6 8 10 12 14 
// R（z）~L1*s+L2*s+L3*s+L4*s+L5*s+L6*s+L7*s 
// （L1到L7的值列在程序中）和
// /|2 14 |-58.45 
// /|
// 注意2s=f-s*f=f-hfsq+s*hfsq，其中hfsq=f*f/2。
// 为了保证log中的错误低于1LP，我们通过
// log（1+f）=f-s*（f-R）（如果f不是太大）
// log（1+f）=f-（hfsq-s*（hfsq+R））计算log。（更准确）
// 
// 3。最后，log（x）=k*Ln2+log（1+f）。
// =k*Ln2_hi+（f-（hfsq-（s*（hfsq+R）+k*Ln2_lo）））
// 这里Ln2被分成两个浮点数：
// Ln2_hi+Ln2_lo，
// 其中n*Ln2_hi总是精确表示| n |<2000。
// 
// 特殊情况：
// 如果x<0（包括-INF），则log（x）为带信号的NaN；
// log（+INF）是+INF；日志（0）是带信号的-INF；
// log（NaN）是没有信号的NaN。
// 
// 精度：
// 根据误差分析，误差始终小于
// 1 ulp（最后一位的单位）。
// 
// 常数：
// 十六进制值是以下
// 常数的预期值。可以使用十进制值，前提是
// 编译器将十进制值转换为二进制值的精度足以产生所示的十六进制值。

// Log返回x的自然对数。
// 
// 特殊情况是：
// Log（+Inf）=+Inf 
// Log（0）=-Inf 
// Log（x<0）=NaN 
// Log（NaN）=NaN 
func Log(x float64) float64 {
	if haveArchLog {
		return archLog(x)
	}
	return log(x)
}

func log(x float64) float64 {
	const (
		Ln2Hi = 6.93147180369123816490e-01 /* 3fe62e42 fee00000 */
		Ln2Lo = 1.90821492927058770002e-10 /* 3dea39ef 35793c76 */
		L1    = 6.666666666666735130e-01   /* 3FE55555 55555593 */
		L2    = 3.999999999940941908e-01   /* 3FD99999 9997FA04 */
		L3    = 2.857142874366239149e-01   /* 3FD24924 94229359 */
		L4    = 2.222219843214978396e-01   /* 3FCC71C5 1D8E78AF */
		L5    = 1.818357216161805012e-01   /* 3FC74664 96CB03DE */
		L6    = 1.531383769920937332e-01   /* 3FC39A09 D078C69F */
		L7    = 1.479819860511658591e-01   /* 3FC2F112 DF3E5244 */
	)

	// 特殊情况
	switch {
	case IsNaN(x) || IsInf(x, 1):
		return x
	case x < 0:
		return NaN()
	case x == 0:
		return Inf(-1)
	}

	// reduce 
	f1, ki := Frexp(x)
	if f1 < Sqrt2/2 {
		f1 *= 2
		ki--
	}
	f := f1 - 1
	k := float64(ki)

	// 
	s := f / (2 + f)
	s2 := s * s
	s4 := s2 * s2
	t1 := s2 * (L1 + s4*(L3+s4*(L5+s4*L7)))
	t2 := s4 * (L2 + s4*(L4+s4*L6))
	R := t1 + t2
	hfsq := 0.5 * f * f
	return k*Ln2Hi - ((hfsq - (s*(hfsq+R) + k*Ln2Lo)) - f)
}
